HiLUK: Scalable Incomplete Factorization Utilizing Combinatorial Methods to Reduce Overheads

Incomplete factorizations are used to approximate the factorization of a sparse coe cient matrix A, such that A = L̄Ū ⇡ LU , and are commonly used as preconditioners for iterative methods, such as GMRES [7]. The approximation is normally achieved by some combination of dropping small value and/or by not allowing fill-in, i.e., zero elements becoming nonzero during factorization, based on levels (k) generated in the elimination tree (ILU-K). Incomplete factorizations are notorious for scaling poorly due to low computational intensity per communication/synchronization (sync). Due to this, very few implementations exist that scale beyond a handful of threads. However, increasing number of light-weight cores require that incomplete factorizations scale in order to not to be the bottleneck in key operations such as preconditioned GMRES. In this work, we present a new incomplete factorization package HiLUK that uses a variety of combinatorial techniques to achieve near linear speedups on x86 and Intel Phi. We only report ILU-0 here for brevity. Sparse factorizations have always required the use of advance combinatorial methods such as graph partitioning and ordering. In order to scale on current systems, a combination of these techniques need to be used to exploit both the matrix sparsity pattern and the underlying hierarchies in modern computer architectures [1]. Synchronization. Traditional methods parallelize incomplete factorization include factoring based on levelsets, coloring, or nested-dissection orderings (ND). However, each of these techniques’ standard implementation requires all threads to sync between levels, colors, or tree levels. These syncs can dominate the execution time as the computational intensity of incomplete factorization is much lower than full factorization. In Figure 1, we present a scatter plot of both the number of rows vs number of syncs required to factor 7 matrices reordered with ND using 16 threads of OpenMP style barriers with level-sets. We see from the plot that the